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We show how to compute the generating function of the self-avoiding polygons on a lattice 
by using the statistical mechanics Schwinger-Dyson equations for the correlation functions of the 
N- vector spin model on that lattice. 
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. I. INTRODUCTION 

It has long been known that nontrivial combinatorial results concerning the enumeration of subgraphs are related 
to the computation of high-temperature (HT) expansions in the statistical mechanics of lattice spin modelsa. 

From a physicist's standpoint there are numerous reasons of interest in the HT expansions. They are powerful 
tools for studying numerically the properties of statistical models and enable us to test directly, in specific cases, the 
validity of the basic assumption of universality in the critical phenomena, as well as to assess the accuracy of the 
estimates of the universal critical observables obtained within the renormalization group theory by such perturbative 
approximations as the e-expansion or the fixed-dimension expansiono. 
t-H , On the other hand, from a mathematician's standpoint, the calculation of long HT expansions poses intriguing 
challenges of combinatorial, algorithmic and analytical nature and, in various cases, it provides generating functions 
£2T" f° r difficult enumeration problems on a latticetl. 

WeJiave devoted this note to a brief discussion of the use of the Schwinger-Dyson (SD) equations in computing HT 
\ seriesBlij for the iV-vector lattice spin modelsE3 and in deriving, as byproducts, combinatorial results concerning the 
■ enumeration of self-avoiding walks(SAW's) on d-dimensional lattices. j-. .— . 

^\ ' The SD equations have been formulated for several lattice spin modelsa O. They are an infinite hierarchy of 
0^ ■ coupled linear equations for the spin correlation functions. We have shown that they can be effectively used as 
recursion relations to compute HT expansions for the correlation functions. In one form or another and under 
different denominations, these or closely related sets of equations have long been known in statistical mechanics and 
in quantum field theory and have been often used in various contexts, in particular recently in lattice field theory (see, 
for example, the very partial list of references^). It seems, however that the SD eqs. have never been systematically 
used to produce HT series, (partial exceptions are Refso), and moreover that they are unfamiliar as a tool of the 



q | combinatorial analysis, so that their potential in this context remains largely unexploited. 

We have formulated and used the SDE eqs. for: i) the A-vccljg ^a^D called 0(A) classical Heisenberg) mode 



on the square and the cubic lattices; ii) the XY (or 0(2)) modeBm3~EJ on the square, on the triangular and on the 
simple cubic lattices; iii) the RP N ~ X model (or generalized Meier-Saupe model)Lj on the square lattice. 

II. THE SD EQUATIONS FOR THE TV- VECTOR MODEL 



The A-vector model is defined as follows: to each site Xi of a lattice (to be definite let us consider a cubic lattice) 
we associate an A-component unit vector s{x*i), called spin, and take as the Hamiltonian 

H{s} = -J2K^)-s^j) (1) 

where (i,j) indicates that the sum extends to nearest neighbor sites only. Therefore 

Z(fi) = f exp{-PH{a})Y[d^Xi) (2) 

J i 

is the partition function and (3 the inverse temperature. [Of course, eq.(^|) and many expressions below have only a 
formal meaning, unless they refer to a finite sublattice A and the limit, as A — > oo, is properly defined. However, the 
rigorous reader can easily work out by himself these well known details that are omitted here to avoid overburdening 
the notation.] If we set 
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(j)(C) = <f>(xi,x 2 ,..x n ;{b itj }) = Y[ ■ s(xj)) b ' J , (3) 

l<i<j<n 

with bi j > 0, the most general correlation function of the model can be written as 

(0(C)) = z{0)J <KC)exp(-l3H{s}) ]J dSfr) (4) 

and therefore it is uniquely identified by the configuration C, namely by the the set of lattice sites x\,x*2i --x n , called 
vertices, and by the integers {bij}, called number of bonds between the vertices (xi,Xj). This notation suggests 
a simple graphical representation for any correlation function (4>(C)) on the lattice. Using the invariance of the 
integration measure in the partition function under a site-dependent O(N) transformation of the spins, we have 
showna how to obtain a convenient set of SD equations for the correlation functions of the model: 

3 n 

(N-2 + 9l mC)) = f3 Yl )) - (HC+))) + (b 12 - 1)<0(C 12 , 12 )) + ]T M#c?ki)} (5) 

^=-3 j=3 

Here g, = ^2j— 2 denotes the number of bonds connected with the vertex x,. The configuration is obtained 
from the configuration C by adding the bond (a?i, x, + e M ), with e M a unit vector in the [i lattice direction, namely 

<t>{C+) = s^) ■ six, + e^4>(C) (6) 

The configuration 0(C~) is obtained from C by changing a bond (afi,^) into (afi + e^,^), namely 



m) = 3{ %t^: { T ] m. (7) 



six,) ■ s(f 2 ) 



The configuration 0(6*12,12) is defined by 



0(Ci2,i 2 ) = ,^^S C }^ s, 2 (8) 

whenever there are at least two bonds connecting the sites x, and X2 in the configuration C. 

Finally, the configuration </>(C^2 y) is obtained from C by removing a bond between ix\,x 2 ), one between ix,,Xj) 
and by adding a new bond between (x 2 ,Xj), namely 

Since C can be any configuration, including the empty configuration C = associated with the trivial correlation 
(0(0)) = 1, eqs. (]^) form an infinite hierarchy of linear (inhomogeneous) equations relating all correlation functions. 

A simple algorithm solves iteratively this system of equations and yields the expansion of any correlation function 
as a power series of the inverse temperature fi. The calculation can be performed on a computer by using exclusively 
integer arithmetics and therefore is completely exact. Moreover it can be formulated in such a way that the coefficients 
of the high temperature expansion are explicitly obtained as ratios of polynomials in the variable N. 

In order to understand the combinatorial relevance of the HT expansion in the N- vector model, it is now sufficient 
to recall the proof£3 that, in the limit as TV — > at fixed (3/N, the correlation function between the spin at the 
origin and the spin at the site x is the generating function of the self-avoiding walks connecting the two sites on the 
lattice, up to a factor g (the coordination number of the lattice). In particular, the problem of the enumeration of 
the self-avoiding polygons with n sides is solved once we determine the (n — 1) — th coefficient in the HT expansion 
of the correlation function between the spins at nearest neighbor sites. The most-extensive results available until 
now for three-dimensional lattices had been obtained by the King's College grouplLEl (up to /3 19 on the simple cubic 
lattice, up to /3 15 on the body-centered cubic lattice and up to /3 13 on the face-centered cubic lattice) and date back 
to more than 25 years ago. We have recently extended these results to order 21 on botL-the simple cubic and the 
body-centered cubic lattices and hope to reach significantly higher orders in a short timetil. Of course, due to the 
enormous progress in computers in the last quarter of century, the traditional direct counting techniques employed by 
the the King's College group would today produce much longer series, but the main point of this note is to present 
an interesting alternative technique rather than an impressive computation. 
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Our extended results are: 

P sc {(3) = 24f3 3 + 264/3 5 + 33 1 2/3 7 + 48 2 40/3 9 + 762096/3 11 + 12673920/? 13 + 

218904768/3 15 + 3891176352/? 17 + 70742410800/3 19 + 1309643747808/3 21 + ... 
on the sc lattice and 

P bcc (f3) = 96/3 3 + 1776/3 5 + 43776/3 7 + 1237920/3 9 + 37903776/3 11 + 1223681760/3 13 + 

41040797376/3 15 + 1416762272736/3 17 + 50027402384640/3 19 + 1799035070369856/3 21 + ... 
on the bcc lattice. 

More results of combinatorial interest emerge, for instance, bwomputing, in the same limit TV — > 0, the suscepti- 
bility, which is the SAW chain generating function etc.. See Ref.lIZl for recent series extensions. 

For the sake of clarity and brevity, we shall not describe here the iterative procedure for the solution of eq. (j^) , but 
ratheij-anly the analogous procedure for the slightly simpler case in which we set N — 2, namely for the so-called XY 
modellla. In this case the site variables are 2-component unit vectors s(xi) — (s^(xi), s^- 2 \xi)) and the Hamiltonian 
and the correlation functions can be parametrized in a much simpler way: 

H{8} = -Y, ■ = - E cos & - 9 i) ( 10 ) 

so that 

Z= / exp(-PH{9})T\d0i (11) 
Jo 

where 6i is the angle formed by the spin at the site £j with a given fixed direction. 
A complete set of correlation functions for this system can be written as follows : 

Ms)) = \ r v( S )exp(-0H{6}) J] Mi (12) 
z Jo 



N 

cp(s) = exp(i 2J QkOk), (13) 
fc=i 

and qk is an integer number that we call the "charge" at the site x^. In this special case of the N- vector model, we 
can put the correlation functions into a one to one correspondence with the "configurations of charges" on the lattice, 
namely 

(ip(s)) <-> s = (xx,q 1 ;x2,q2; ...;x N ,q N ). 

Let us observe that the Hamiltonian has the global symmetry 9i — > 9i + tp. On two-dimensional lattices, by the 
Mcrmin- Wagner theorem, the continuous symmetries cannot be broken, therefore 

N 

Ms))^0 5> fc = 0, (14) 

k=l 

at any temperature. This means that nontrivial correlation functions are associated only to the globally "neutral" 
charge configurations. On three-dimensional lattices, this is still true at high temperatures, while the symmetry is 
broken at low temperatures. 

For this model the simplest set of SD eqs. is the infinite hierarchy 

/ U l d8 l — [^(s)exp(-f3H{9})}=0 (15) 
Jo d9 k 

for all possible <f(s) and taking Xk G s. These identities are an obvious consequence of the 2tt periodicity of the 
quantity in square parentheses. More explicitly we have 
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Ms)) = f- ((exp[-i(e h - 6 k+fi )]<p(s)) - {exp[i{9 k - 6 k+ »)]<p{s))) (16) 

where the sum extends to all of the nearest neighbor sites of Xk ■ Introducing the notation s+ (s~ ) for the lattice 
charge configuration obtained increasing (decreasing) by 1 the charge at the site x k and decreasing (increasing) by 1 
the charge at the site x k + e^, we can write eqs. ( |16| ) more compactly as follows: 

(v( 8 )> = ^E(<^)>-m<)>)- ( 17 ) 

For example, on the square lattice let us consider the charge configuration s = (0, 0, 1; 1, 0,-1), so that (tp(s)) is 
the elementary nearest neighbor correlation function (exp[i(9oo — #io)])- Let us also assume that the derivative 
acts on the first site of each lattice charge configuration, namely take k = 1, then s+ and s~ represent the following 
four configurations 

s+ ee {(0,0, 2; 1,0, -2), (0,0, 2; 1,0,-1; 0,1,-1), (0,0, 2; 1,0,-1; -1,0,-1), 

(0,0,2;1,0,-1;0,-1,-1)} (18) 



b- = {0,(0, 1,1; 1,0,-1), (-1,0,1; 1,0,-1), (0,-1,1; 1,0,-1)} (19) 

It is easy now to understand intuitively how to produce a HT expansion: we have only to notice that in the previous 
example the correlation function (exp[i(8oo — 0io)}) on the lhs of the SD eq. has a HT expansion beginning at order 
/3, while on the rhs there is a linear combination of correlation functions whose HT expansions begin at order /J 2 , 
except for (tp(0)) — 1. Therefore we obtain 

(exp{t(e o-e w )}} = ^(l + 0({3 2 )). (20) 

This remark immediately yields the first term of the HT expansion for the nearest neighbor correlation function. The 
successive term, (which is actually 0(/3 3 )), is obtained by the same mechanism using the SD eqs. for the correlation 
functions (<p(s~)) and (ip(s+)) in the O(0 2 ) terms on the rhs and so on. Once this simple mechanism has been made 
clear, it is convenient to formulate the procedure in a way which also suggests how to make it completely automatic 
and therefore ready to be translated into a computer code. 
To this purpose let us first define a linear space of states as follows: 

Consider the set of all possible distinct neutral charge configurations on the lattice. Two such charge configurations 
s ee (a?i, gx! %2i Q2\ • %N, Qn) an d s' ee (x[, q[; q^; ...; x* M , <^ M ) are considered distinct if they cannot be identified 
by a permutation of labels and/or by a lattice roto-translation and/or by changing the sign of all charges. To each 
charge configuration s, we associate a state vector represented by a ket |s); to the empty lattice s = 0, we associate 
the ket |0). We define a scalar product as follows: (s|s') = 1 if s = s' and (s|s') = otherwise. The set of all possible 
states is therefore an orthonormal basis of an infinite-dimensional linear vector space on which we can define a linear 
operator K as follows: 

#10) =0 (21) 

*I-> = ^£(M-K!:>) C 22 ) 

here q\ is the charge at the site x\. 
Now consider the identity operator 

1 = (1 - PK)(1 - (3K)- 1 (23) 
and take its matrix element between the states (0| and |s), namely: 

(0|(l-^)(l-/3#)- 1 |s)=0. (24) 

From ( p4| ) it follows : 
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(0|(1 - PKY 1 ^) = /9<0|(1 - PKV'K |s>. (25) 
Using the definition of X (|2l|),(p2[) we have finally: 



(OKl-^)- 1 ^) = ^;I>|(1 -^JO-^ls-) - |s+)). (26) 

We conclude that (0|(1 — /?A') _1 |s) satisfies the same SD eq. as (y(s)) and therefore 

(^(s)) = {0|(l-/5X)- 1 |s). (27) 

Now the solution of the SD eqs. can simply be obtained by a power series expansion of the resolvent operator 
(1-PK)- 1 . We have 

Ms))=£/3»(0|X"|s>. (28) 
n 

This is precisely the HT expansion of ((p(s)). The n-th order coefficient 

c n = (0\K n \s)= £ <0|A'|s 1 }<s 1 |A> 2 }...<s n _ 1 |A» (29) 

is a particular matrix element of the n-th power of the matrix (s|A"|s'}. 

The structure of the matrix K is quite simple: it is a (very) sparse matrix with elements 

(sills') = ^ (30) 
2qi 

if s = s'^ and otherwise. 

Therefore in every column there appear at most 2g nonvanishing elements (g being the coordination number of the 
lattice). 

In order to calculate a HT series up to some maximal order j3 M , it is sufficient to compute the matrix A' on a finite 
truncated basis. We can then implement the following iterative procedure for computing the HTE of the correlation 
function ((f(s)) up to the order M. 

In order to be definite we shall continue to consider the case of a square lattice. 

First iteration: 

consider the state |s) and call it state Compute K|l). This produces the 8 new states |sj) = (|2), |3), .., |9)). 
Memorize the matrix elements of column 1, namely: 
K 21 = {2\K\1),K 31 = {3\K\1),...,K ?1 = {9\K\1). 

Clearly these are the only nonvanishing elements of the first column of the matrix K. 

(Remember that the operator K acts only on the "first" charge of each state, which can be arbitrarily chosen.) 
Second iteration: 

Act with K on the new states |2), |3), .., |9). Then 8 2 new configurations are produced: 10) , |11), .., |73). Memorize 
now the nonvanishing elements of columns 2,3, ..9 namely ATio,2, ^11.2, ^73,9- 
Continue until M iterations are completed. 

This is a correct, but definitely naive and inefficient procedure: the growth of the truncated basis on which we 
compute the elements of the matrix AT is very fast. It appears that, in order to perform the (M+l)-th iteration, we 
should store 8 M states ( or (2g) for a general lattice with coordination number g). Therefore M ~ 10, or even less 
than that, would seem to be an insuperable upper limit for this kind of computation. However, only a small subset of 
the states produced in the iteration are distinct with respect to the natural equivalence relation mentioned above. We 
must obviously exploit this fact by always identifying equivalent states, namely by working with equivalence classes 
of states. This reduces drastically the growth of the truncated basis. For example, on the square lattice, the growth is 
reduced from 8 M to some (~ 3) M . Another very useful observation is that, after the [M/2]-th iteration, an increasing 
fraction of the new states produced do not actually contribute to the matrix elements of (0|AT r |s), provided that 
r < M. Therefore these states are irrelevant for the calculation and can be simply discarded. This entails a further 
reduction of the dimension of the truncated basis to some more manageable (~ 2) M . These simple remarks show why 
much higher orders than M ~ 10 can be attained. 

We have sketched the simplest recursive approach to the SD eqs. (||). The procedure needed for the general 
iV-vector model is only slightly more intricate and the interested reader can find all details in Ref.B. 
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III. CONCLUSIONS 



In conclusion, we stress that conceptual simplicity is the main advantage of the SD method. It enables us to 
bypass the usual elaborate graphical formulation of the HT expansion procedures in statistical mechanics, namely the 
construction and selection of the relevant graphs, the computation of the graph embedding constants in the lattice 
etc., as well as any other laborious combinatorial inputs. Indeed all correlation functions of the system under study 
are produced directly by an elementary recursion. Moreover the method is very flexible and widely applicable. This 
simplicity also entails the limitations of the technique: after iterating sufficiently many times the recursion relations, 
a plethora of contributions is generated and, for each given model, it takes some ingenuity to classify them and to 
organize efficiently the computation by managing carefully and economically the storage resources of the computer. 
Moreover, due to presence of a summation over \i in eqs. @) and (0), the number of contributions at a given 
iteration increases rapidly with the dimensionality of the lattice. In our first tests of the method, the maximal order 
of expansion was limited only by the size of the computer memory, (which, only a few years ago, was usually two orders 
of magnitude smaller than it is today in a common pc), whereas the required CPU times were always very modest. For 
this reason, after performing a few computations, we had set aside the recursion methods and turned to the traditional 
linked-cluster-expansion techniquesLD. We believe however, that the Schwinger-Dyson technique will be increasingly 
useful and easy to use as larger mass memories become available in computers. Certainly, the very impressive results 
obtained by algebraic techniques^ in low-dimensional models with discrete site- variables (for instance SAW and Ising 
models) are unlikely to be improved or even to be reproduced by these procedures. The recursion approach, however, 
remains superior for models with continuous variables, like the XY model, where the algebraic methods become 
inefficient. Therefore, we believe that the Schwinger-Dyson method deserves further exploration in a wider context, 
both from the formal and from the algorithmic point of view. 
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